package filtros;

/***********************************************************************************
FiltroShenCastan.java - Implementa o filtro de Shen e Castan para detec��o de 
												bordas.

Direitos Autorais Reservados (c) 2008 Embrapa Inform�tica Agropecu�ria

Este programa � software livre; voc� pode redistribu�-lo e/ou modific�-lo sob os
termos da Licen�a P�blica Geral GNU conforme publicada pela Free Software
Foundation; tanto a vers�o 2 da Licen�a, como (a seu crit�rio) qualquer vers�o
posterior.

Este programa � distribu�do na expectativa de que seja �til, por�m, SEM NENHUMA
GARANTIA; nem mesmo a garantia impl�cita de COMERCIABILIDADE OU ADEQUA��O A UMA
FINALIDADE ESPEC�FICA. Consulte a Licen�a P�blica Geral do GNU para mais detalhes.

Voc� deve ter recebido uma c�pia da Licen�a P�blica Geral do GNU junto com este
programa; se n�o, escreva para a Free Software Foundation, Inc., no endere�o 59
Temple Street, Suite 330, Boston, MA 02111-1307 USA.

**********************************************************************************

Parte integrante do Comunicado T�cnico 089/08.
Acesso: <http://www.cnptia.embrapa.br/files/comtec89.pdf>

Descri��o:

	Implementa o filtro de detec��o de bordas de Shen e Castan.
 
	A implementa��o do filtro tem por refer�ncia o artigo:

	SHEN, J.; CASTAN, S. An optimal linear operator for step edge detection. 
	   Computer Vision, Graphics, and Image Processing: Graphical Models and 
		 Understanding. 54(2):112-133, 1992.

Par�metros de entrada:

  1. imagem			=>	nome do arquivo imagem.
	2. fatorB (b)	=>	fator de suaviza��o
	3. windowSize	=>	largura da janela do filtro (kernel)

Nota: Para um resultado razo�vel, usar:
			fatorB = 0.75 (fator de suaviza��o 0 < b < 1; equa��o (3)) 
			windowSize = 11

Sa�da:

  Painel com a imagem original e imagem com as bordas.

Desenvolvido por:

    Jos� Iguelmar Miranda & Jo�o Camargo Neto

Informa��es do CVS:
       $Source$:
       $Revision$:
       $Date$:

***********************************************************************************/

// pacotes gen�ricos
import java.io.File;

// pacotes do AWT
import java.awt.image.WritableRaster;
import java.awt.image.BufferedImage;
import java.awt.image.Raster;
import java.awt.GridLayout;

// pacote Swing para interface gr�fica
import javax.swing.JFrame;
import javax.swing.JLabel;
import javax.swing.ImageIcon;
import javax.swing.JScrollPane;

// leitura em modo imediato - para J2SE 1.4+
import javax.imageio.ImageIO;

public class FiltroShenCastan extends JFrame {

	// Atributos
	private static int OUTLINE = 25;

  public static void main(String args[]) {

		// Par�metros do operador ShenCastan
		double b = 0.8;	
		int windowSize = 7;

               /* args[0] = "C:\Users\Leandro\Desktop\testeImagem.jpeg";
                args[1] = "0.75";
                args[2] = "11";*/
		// Verifica passagem de par�metros
    if (args.length != 3) {
			String msga = "\nUso: java -cp . FiltroShenCastan <imagem>";
			String msgb = " <fatorB><windowSize>";
			System.out.println(msga + msgb);
      System.exit(0);
    }

		// Mostra JFrame decorado pelo Swing
		JFrame.setDefaultLookAndFeelDecorated(true);

		try	{
	    b = Double.parseDouble(args[1]);
		} catch (NumberFormatException e)	{
			System.out.println("Erro no valor de 'fatorB'");
      System.exit(0);
		}

		// b � um fator de suaviza��o
	  if (b < 0)
			b = 0.0;
    else 
			if (b > 1.0)
				b = 1.0;

		try	{
	    windowSize = Integer.parseInt(args[2]);
		} catch (NumberFormatException e)	{
			System.out.println("Erro no valor de 'windowSize'");
      System.exit(0);
		}

		System.out.println("\nFiltro ISEF - Parametros:\n");
		String st = "\n\tFator suavizacao: ";
	  System.out.println(st + b);
		st = "\n\tLargura da janela para gradiente adaptativo: ";
	  System.out.println(st + windowSize);

		// Tempo inicial
		long eq_time = System.currentTimeMillis();

    // Chama o construtor
		FiltroShenCastan s = new FiltroShenCastan(args[0], windowSize, b); 

		// Tempo final
    eq_time = System.currentTimeMillis() - eq_time;
		String msg = "\n\tShen-Castan: tempo de execucao ";
    System.out.println(msg + eq_time + " milisseg.");

		// Encerra a aplica��o clicando no "close"
		s.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
		s.setVisible(true);
	}
  
	// Construtor
	public FiltroShenCastan(String aFile, int windowSize,	double b) {

		// Define objeto BufferedImage para encapsular a imagem
    BufferedImage im = null, molduraImagem = null, 
									imagemSemMoldura = null, res = null;
		int w, h;
		JLabel img1, img2;
      
    // Certifica que a imagem existe
    File file = new File(aFile);
    try {
        im = ImageIO.read(file);
    } catch(Exception e) {
        System.out.println("Imagem '" + aFile + "' nao existe.");
        System.exit(0);
    }

    // Atribui nome ao frame
    setTitle("Shen-Castan: " + file.getName());

		// Processamento do filtro
		molduraImagem = funcaoMoldura(im, OUTLINE);

    WritableRaster resWR = molduraImagem.copyData(null); 
	  res = new BufferedImage(molduraImagem.getColorModel(), resWR, false, null);
		res = shenOp(molduraImagem, res, windowSize, b);
		imagemSemMoldura = funcaoSemMoldura(res, OUTLINE);

		w = im.getWidth();
		h = im.getHeight();

		// Define GUI com objetos do Swing
    getContentPane().setLayout(new GridLayout(1, 2));
    img1 = new JLabel(new ImageIcon(im));
    img2 = new JLabel(new ImageIcon(imagemSemMoldura));
		getContentPane().add(new JScrollPane(img1));
		getContentPane().add(new JScrollPane(img2));
    setSize(2*w, h);
	}

	// Retorna uma BufferedImage envolta em uma moldura de largura width
	private BufferedImage funcaoMoldura(BufferedImage im, int width) {
		int iP, jP, tipo, colunas, linhas, col, lin, pixel;
		BufferedImage newImg;

		colunas = im.getWidth();
		linhas = im.getHeight();
		width += 2;
    tipo = BufferedImage.TYPE_BYTE_GRAY;
    newImg = new BufferedImage(colunas + 2*width, linhas + 2*width, tipo);
    WritableRaster newImgWR = newImg.getRaster();
    Raster imR = im.getRaster(); 

		col = newImg.getWidth();
		lin = newImg.getHeight();

		for (int i = 0; i < lin; i++)
		  for (int j = 0; j < col; j++)
			  newImgWR.setSample(j, i, 0, 0);

		for (int i = width; i < lin - width; i++)
		  for (int j = width; j < col - width; j++) {
		    jP = j - width;
		    iP = i - width;
				pixel = imR.getSample(jP, iP, 0);
			  newImgWR.setSample(j, i, 0, pixel);
			}

		return newImg;
	}

	// M�todo para realizar opera��es do filtro
	private BufferedImage shenOp(BufferedImage im, BufferedImage res,
															 int windowSize, double b) {

		double[][] imagemFloat, smoothedImg;
		int colunas, linhas;
		int[][] imagemInt = null;
		BufferedImage imgBinLaplace = null;		// imagem bin�ria do Laplaciano
   
    Raster imR = im.getRaster(); 
    Raster resR = res.getRaster(); 
    WritableRaster resWR = res.getRaster();

		colunas = im.getWidth();
		linhas = im.getHeight();
		imagemInt = new int[colunas][linhas];

		// Converte a imagem de entrada para ponto flutuante
		imagemFloat = new double[colunas][linhas];
		for (int i = 0; i < linhas; i++)
			for (int j = 0; j < colunas; j++)
				imagemFloat[j][i] = (double)imR.getSample(j, i, 0);

		// Suaviza a imagem de entrada usando, recursivamente, o filtro
		// ISEF
		smoothedImg = computeISEF(imagemFloat, b);

		// Calcula imagem bin�ria do Laplaciano (ibl) a partir da imagem suavizada.
		// imgBinLaplace � uma imagem bin�ria {0, 1}.
		imgBinLaplace = computeIBL(smoothedImg, imagemFloat); 

		// Faz detec��o de bordas usando ibl com limiar de gradiente
		// Esta fun��o altera apenas os valores de 'imagemFloat', que 
		// armazena os valores das bordas (gradiente).
		imagemFloat = locateZeroCrossings(smoothedImg, imgBinLaplace, windowSize);

		// Limiariza bordas
		for (int i = 0; i < linhas; i++)
		  for (int j = 0; j < colunas; j++)
			  if (imagemFloat[j][i] > 0)
					imagemFloat[j][i] = 0;
				 else
					imagemFloat[j][i] = 255;

		// Converte o imagemFloat double para int.
		for (int i = 0; i < linhas; i++)
		  for (int j = 0; j < colunas; j++)
				imagemInt[j][i] = (int)imagemFloat[j][i];

		// Cria imagem com o imagemInt de inteiros...
		res = criaImagem(imagemInt);

		// ... e retorna a imagem.
		return res;
	}
	
	private BufferedImage criaImagem(int[][] in) {
		int w, h, tipo, pixels[], ind = 0;
		BufferedImage dest = null;

	  w = in.length;
		h = in[0].length;
		pixels = new int[w*h];
    tipo = BufferedImage.TYPE_BYTE_GRAY;
    dest = new BufferedImage(w, h, tipo);
		WritableRaster destWR = dest.getRaster();

		for(int lin = 0; lin < h; lin++)
			for(int col = 0; col < w; col++)
				pixels[ind++] = in[col][lin];

		destWR.setPixels(0, 0, w, h, pixels);
		return dest;		
	}

	// Filtro recursivo de ISEF
	// Shen & Castan, CVIGP March 1992
	// x = imagem no formato de ponto flutuante, conforme equa��o (3).
	// y = imagem resultado, suavizada.
	private double[][] computeISEF(double[][] imgFloat, double b) {

		int nLinhas, nColunas;
		double[][] causalComp, antiCausalComp, y;
   
		nColunas = imgFloat.length;
		nLinhas = imgFloat[0].length;

		// Armazena componente causal
		causalComp = new double[nColunas][nLinhas]; 
		
		// Armazena componente anti-causal
		antiCausalComp = new double[nColunas][nLinhas]; 
   
		// Primeiro, aplica o filtro na dire��o vertical (nas linhas).
		// O resultado do filtro � armazenado na vari�vel y.
		y = aplicaISEFVertical(imgFloat, causalComp, antiCausalComp, b);

		// Segundo, aplica o filtro na dire��o horizontal (nas colunas),
		// usando a imagem filtrada nas linhas, y.
		y = aplicaISEFHorizontal(y, causalComp, antiCausalComp, b); 
		return y;
	}

	// Filtra linhas da imagem de entrada, x, e coloca resultado em y.
	private double[][] aplicaISEFVertical(double[][] x, double[][] causalC, 
											 double[][] antiCausalC, double b) {

		double b1, b2, y[][];
		int nLinhas, nColunas;

		nColunas = x.length;
		nLinhas = x[0].length;
		y = new double[nColunas][nLinhas];
   
		b1 = (1.0 - b)/(1.0 + b);
		b2 = b*b1;
   
		// Calcula as condi��es de bordas
		for (int col = 0; col < nColunas; col++) {
			// Borda existe para primeira e �ltima colunas
			causalC[col][0] = b1 * x[col][0];	
			antiCausalC[col][nLinhas-1] = b2 * x[col][nLinhas-1];
		}

		// Calcula componente causal
		for (int lin = 1; lin < nLinhas; lin++)
			for (int col = 0; col < nColunas; col++)
				causalC[col][lin] = b1 * x[col][lin] + b * causalC[col][lin-1];

		// Calcula componente anti-causal
		for (int lin = nLinhas - 2; lin >= 0; lin--)
			for (int col = 0; col < nColunas; col++)
				antiCausalC[col][lin] = b2 * x[col][lin] + b * antiCausalC[col][lin+1];

		// Caso de borda para calcular sa�da do primeiro filtro
		for (int col = 0; col < nColunas-1; col++)
			y[col][nLinhas-1] = causalC[col][nLinhas-1]; 

		// Calcula a sa�da do primeiro filtro e armazena em y
		// Equivale � soma dos componentes causal e anti-causal
		for (int lin = 0; lin < nLinhas - 2; lin++)
			for (int col = 0; col < nColunas - 1; col++)
				y[col][lin] = causalC[col][lin] + antiCausalC[col][lin+1];
		return y;
	}  

	// Filtra colunas da imagem de entrada, x, e coloca resultado em y.
	private double[][] aplicaISEFHorizontal(double[][] x, double[][] causalC, 
												 double[][] antiCausalC, double b) {
	
		double b1, b2, y[][];
		int nLinhas, nColunas;

		nColunas = x.length;
		nLinhas = x[0].length;
		y = new double[nColunas][nLinhas];
     
		b1 = (1.0 - b)/(1.0 + b);
		b2 = b*b1;
   
		// Calcula as condi��es de bordas
		for (int lin = 0; lin < nLinhas; lin++) {
			causalC[0][lin] = b1 * x[0][lin];
			antiCausalC[nColunas-1][lin] = b2 * x[nColunas-1][lin];
		}

		// Calcula componente causal
		for (int col = 1; col < nColunas; col++)
			for (int lin = 0; lin < nLinhas; lin++)
				causalC[col][lin] = b1 * x[col][lin] + b * causalC[col-1][lin];

		// Calcula componente anti-causal
		for (int col = nColunas - 2; col >= 0; col--)
			for (int lin = 0; lin < nLinhas; lin++)
				antiCausalC[col][lin] = b2 * x[col][lin] + b * antiCausalC[col+1][lin];

		// Caso de borda para calcular sa�da do primeiro filtro
		for (int lin = 0; lin < nLinhas; lin++)
			y[nColunas-1][lin] = causalC[nColunas-1][lin];

		// Calcula a sa�da do segundo filtro e armazena em y
		// Equivale � soma dos componentes causal e anti-causal
		for (int lin = 0; lin < nLinhas; lin++)
			for (int col = 0; col < nColunas - 1; col++)
				y[col][lin] = causalC[col][lin] + antiCausalC[col+1][lin];
		return y;
	}  

	// Calcula o Laplaciano limitado por banda da imagem de entrada
	private BufferedImage computeIBL(double[][] smoothed, double[][] imgFloat) {

		int tipo, pixel, nLinhas, nColunas;
		BufferedImage imgBinLaplace = null;
    
		nColunas = smoothed.length;
		nLinhas = smoothed[0].length;

		tipo = BufferedImage.TYPE_BYTE_GRAY;
    imgBinLaplace = new BufferedImage(nColunas, nLinhas, tipo);
    WritableRaster bliWR = imgBinLaplace.getRaster();

		for (int lin = 0; lin < nLinhas; lin++)
			for (int col = 0; col < nColunas; col++)
				bliWR.setSample(col, lin, 0, 0);
  
		// A ibl � calculada tomando a diferen�a entre a imagem suavizada
		// e a imagem original (imgFloat), conforme equa��o (4) do texto do CT.
		// A imagem ibl � criada atribuindo todos os valores em ibl com 1 onde o 
		// Laplaciano � positivo e 0, em caso contr�rio.
		for (int lin = 0; lin < nLinhas; lin++)
			for (int col = 0; col < nColunas; col++) {
				// Ignora os pixels ao redor das bordas da imagem.
				if (lin < OUTLINE || lin >= nLinhas - OUTLINE ||
            col < OUTLINE || col >= nColunas - OUTLINE)
					continue;
				if ((smoothed[col][lin] - imgFloat[col][lin]) > 0.0)
					pixel = 1;
				else
					pixel = 0;
				bliWR.setSample(col, lin, 0, pixel);
	  }

		return imgBinLaplace;
	}

	// Esta fun��o muda os valores de 'orig.'
//	private void 
//		locateZeroCrossings(double[][] orig, double[][] smoothed,
//												BufferedImage bli, int windowSize) {
	private double[][] 
		locateZeroCrossings(double[][] smoothed, BufferedImage bli, int windowSize) {

		int nLinhas, nColunas;
		double grad;
		double[][] orig;

		nColunas = smoothed.length;
		nLinhas = smoothed[0].length;

		orig = new double[nColunas][nLinhas];

		for (int lin = 0; lin < nLinhas; lin++) {
			for (int col = 0; col < nColunas; col++) {
				// Ignora os pixels ao redor das bordas da imagem.
				if (lin < OUTLINE || lin >= nLinhas - OUTLINE ||
						col < OUTLINE || col >= nColunas - OUTLINE) {
					orig[col][lin] = 0.0;
				} else
					// Em seguida, verifica se o pixel � um cruzamento zero do Laplaciano
					if (cruzamentoZero(bli, smoothed, lin, col)) {
						// Realiza limiariza��o com gradiente
						grad = calculaGradAdaptativo(bli, smoothed, lin, col,
										windowSize);
						orig[col][lin] = grad;
					} else
						orig[col][lin] = 0.0;		    
			}
		}
		return orig;
	}

	// Acha os cruzamentos zero do Laplaciano (buff).
	// orig = imagem suavizada.
	private boolean 
		cruzamentoZero(BufferedImage buff, double[][] orig,	int lin, int col) {

		// Testa se � cruzamento zero (cz) do Laplaciano. Certifica se o princ�pio 
		// da correspond�ncia do sinal do cz � satisfeito:
		// um sinal positivo do cz deve ter a primeira derivada positiva. Isso 
		// significa que a segunda derivada muda de positivo para negativo � 
		// medida que se passa pela borda.
    Raster buffR = buff.getRaster();

		// cz positivo
		if (buffR.getSample(col, lin, 0) == 1 && 
				buffR.getSample(col, lin + 1, 0) == 0) { 
			if (orig[col][lin+1] - orig[col][lin-1] > 0)
				return true;
			else
				return false;
		} else
			// cz positivo
			if (buffR.getSample(col, lin, 0) == 1 && 
					buffR.getSample(col+1, lin, 0) == 0 ) {
				if (orig[col+1][lin] - orig[col-1][lin] > 0)
					return true;
				else 
					return false;
			} else
				/* negative z-c */
				if (buffR.getSample(col, lin, 0) == 1 && 
						buffR.getSample(col, lin - 1, 0) == 0) {
					if (orig[col][lin+1] - orig[col][lin-1] < 0)
						return true;
					else
						return false;
				} else
					// cz negativo
					if (buffR.getSample(col, lin, 0) == 1 &&
							buffR.getSample(col-1, lin, 0) == 0 ) {
						if (orig[col+1][lin] - orig[col-1][lin] < 0)
							return true;
						else
							return false;
					} else
						// n�o � cz
						return false;
	}

	private double 
		calculaGradAdaptativo(BufferedImage BLI_buffer, double[][] orig_buffer, 
													int lin, int col, int windowSize) {

		double somaBorda, somaFora, mediaBorda, mediaFora;
		int quantBorda, quantFora;
   
		somaBorda = somaFora = 0.0;
		quantBorda = quantFora = 0;
    Raster buffR = BLI_buffer.getRaster();

		for (int i = (-windowSize/2); i <= (windowSize/2); i++) {
			for (int j = (-windowSize/2); j <= (windowSize/2); j++) {
				if (buffR.getSample(col+j, lin+i, 0) == 1) {
	        somaBorda += orig_buffer[col+j][lin+i];
	        quantBorda++;
	     } else {
	        somaFora += orig_buffer[col+j][lin+i];
	        quantFora++;
	     }
	   }
		}
		if (somaFora != 0.0)
			mediaFora = somaFora/quantFora;
		else
			mediaFora = 0.0;
		if (somaBorda != 0.0) 
			mediaBorda = somaBorda/quantBorda;
		else
			mediaBorda = 0;
   
		return (mediaFora - mediaBorda);
	}

	private BufferedImage funcaoSemMoldura(BufferedImage im, int width) {

		int w, h, nr, nc, tipo;
		BufferedImage old = null;

		w = im.getWidth();
		h = im.getHeight();
		width += 2;

    tipo = BufferedImage.TYPE_BYTE_GRAY;
		old = new BufferedImage(w - 2*width, h - 2*width, tipo);
		WritableRaster oldWR = old.getRaster();
		Raster imR = im.getRaster();

		nc = old.getWidth();
		nr = old.getHeight();
		for (int i = 0; i < nr - 1; i++) {
			for (int j = 1; j < nc; j++) {
				oldWR.setSample(j, i, 0, imR.getSample(j+width, i+width, 0));
				oldWR.setSample(j, nr-1, 0, 255);
			}
			oldWR.setSample(0, i, 0, 255);
		}

		return old;
	}

}
